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ABSTRACT 

Differential rotation between the neutron star crust and a more rapidly rotating 
interior superffuid leads to frictional heating that affects the star's long-term thermal 
evolution and resulting surface emission. The frictional heating rate is determined 
by the mobility of the vortex lines that thread the rotating superfluid and pin to the 
inner crust lattice. If vortex pinning is relatively strong, a large velocity difference 
develops between the inner crust superfluid and the crust, leading to a high rate 
of heat generation by friction. Here we present the results of thermal evolution 
simulations based on two models of the vortex pinning forces that bracket a range of 
plausible pinning strengths. We include the effects of superfluidity, magnetic fields, 
and temperature gradients. As representative standard and accelerated neutrino 
emission processes taking place in the core, we consider the modified Urea process in 
normal baryonic matter, and the much faster quark Urea process. Comparison of our 
results with neutron star surface temperature data, including the recent temperature 
measurement of the Geminga pulsar, shows that stars with soft equations of state 
and modest frictional heating are in closest agreement with the data; stars with 
stronger frictional heating have temperatures inconsistent with the upper limit of PSR 
1929+10. Stiffer stars undergoing standard cooling generally have temperatures lying 
above the Vela detection, a situation worsened by the inclusion of frictional heating. 
Stars undergoing accelerated cooling without frictional heating have temperatures 
that fall far below most temperature measurements; the Vela and Geminga detections 
being the most compelling examples. Only in stiff stars, which have thick crusts, can 
the inclusion of strong frictional heating raise the temperature at late stages in the 
evolution to a level consistent with the data. However, such a large amount of heating 
leads to a temperature at ~ 1000 yr in excess of the Crab upper limit. Suppression of 
accelerated neutrino emission processes, perhaps by superfluid pairing in the core, may 
yield acceptable cooling models. 
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1. INTRODUCTION 

The study of neutron star thermal evolution affords a means of probing the properties of dense 
matter. A newborn neutron star is very hot, with a temperature well above an MeV. The star cools 
initially through neutrino emission, and later (after ~ 10 4 yr) cools primarily through the emission 
of photons from the surface. The neutrino emission process occurring in the core depends on the 
state of matter there. In a core composed of normal baryons, neutrino emission can occur through 
the relatively slow modified Urea process ("standard" cooling), or, if the core proton fraction is 
high enough, through the much faster direct Urea process ("accelerated" cooling). Possible states 
of the core include condensed pions or kaons, bulk quark matter, or hyperon-rich matter; neutrino 
emission processes in each of these possible states would have accelerated rates comparable to the 
direct Urea process. Superfluid pairing of the nucleons is expected to occur below temperatures 
of ~ 1 MeV, and could have a dramatic effect on neutron star thermal evolution by reducing the 
matter's specific heat, and, depending on the core emission processes that actually occur, the 
neutrino emissivity. Furthermore, the presence of a neutron superfluid in the inner crust would 
lead to frictional heating, possibly raising the star's surface temperature significantly. 

Generally, stellar evolution models simulating standard cooling (i.e., dominated by the 
modified Urea process) based on moderately stiff to very stiff equations of state predict surface 
temperatures greater than observed; whereas models simulating the accelerated neutrino emission 
from an exotic or proton-rich core predict late-time surface temperatures that are far too low 
(see e.g., Friman & Maxwell 1979 and Van Riper 1991). Internal heating, associated with friction 
between the normal matter of the inner crust and a more rapidly rotating neutron superfluid, 
might allow agreement with data for models involving accelerated cooling, while exacerbating the 
problem in models of standard cooling. 

Frictional heating has been examined by several authors. Alpar et at. (1987) assumed an 
approximate balance between frictional heating and photon emission from the stellar surface, and 
compared the predicted surface temperature with the EXOSAT upper limit on the temperature for 
PSR 1929+10; they concluded that only a small degree of frictional heating is allowed. Shibazaki 
&; Lamb (1989) explored the effects of internal heating on the thermal evolution of isothermal stars 
older than 10 3 yr, including possible magnetic field evolution. Taking the internal heating rate as 
a free parameter, they compared their results with EINSTEIN and EXOSAT detections and upper 
limits. They found standard cooling in a stiff star to be inconsistent with the Vela detection, 
and consistent with other pulsars for various values of the heating rate. Accelerated cooling was 
found to predict temperatures falling below the surface temperatures of several pulsars, even for 
unrealistically high heating rates. Van Riper (1991a) considered constant internal energy sources, 
as would occur for steady state frictional heating, in nonisothermal neutron stars. Van Riper, 
Epstein, & Miller (1991) showed that the heat generated during a pulsar glitch could raise the 
stellar surface temperature over time scales of hours to months. Umeda, Shibazaki, Nomoto, 
& Tsuruta et al. (1993, hereafter USNT) used the parametrized heating rate of Shibazaki & 
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Lamb (1989) in nonisothermal simulations, and compared their results with ROSAT data that 
has recently become available. They found that standard cooling in stiff to medium stars gives 
temperatures that exceed those of PSR 0656+14 and Vela. Moreover, strong heating in a stiff 
star was ruled out by PSR 1055-52, and possibly by the PSR 1929+10 upper limit. The study of 
USNT is the most complete investigation to date, and we compare our results to theirs below. 

The purpose of our study is to consider the consequences of inner crust frictional heating on 
neutron star thermal evolution by using a microscopic description of the superfluid coupling to the 
crust, rather than a parametrized form. In our analysis we use the most recent ROSAT data of 
cooling neutron stars, including the Geminga pulsar detection. We investigate a standard cooling 
model in which the modified Urea process on nucleons is the dominant cooling process in the core. 
As an example of an accelerated emission mechanism, we consider the quark Urea process. Much 
recent work has focused on the microscopic physics determining the coupling between the inner 
crust superfluid and the rest of the star (Alpar 1977; Pines et al. 1980; Alpar et al. 1984; Epstein 
& Baym 1988; Bildsten & Epstein 1989; Link & Epstein 1991, hereafter LE; Link, Epstein, & 
Baym 1993, hereafter LEB; Chau & Cheng 1993a, 1993b). In this study, we examine internal 
heating using the microscopic description for the neutron superfluid dynamics developed in LE 
and LEB. In this theory, the frictional heating rate depends on the stellar temperature and the 
strength of the interactions that pin superfluid vortex lines to the stellar crust. 

Our simulations show that accelerated cooling with frictional heating is inconsistent with 
the data for both young and old pulsars; internal heating occurring at a rate consistent with the 
temperatures of old pulsars produces overheating in the early states of evolution. Standard cooling 
with modest frictional heating of a very soft star (i.e., with a maximum mass of near 1.4 M Q ) is 
marginally consistent with all current temperature data. Stronger frictional heating produces a 
temperature at the age of PSR 1929+10 in excess of the upper limit for that pulsar. 

This paper is organized as follows. In §2 we discuss the input physics of our models; we derive 
the frictional heating rate, and describe the two vortex pinning models used. In §3 we describe our 
numerical simulations. We compare our results with those of other authors and with data from the 
EINSTEIN, EXOSAT and ROSAT satellites in §4, and in §5 we present our conclusions. Appendix 
A gives the vortex creep equations from the theory of LEB. Appendix B contains a comparison of 
our time-dependent source calculations with simulations of cooling with time-independent sources. 



2. INPUT PHYSICS 

2.1. Cooling Processes 

A neutron star cools by neutrino emission from the interior and photon emission from the 
surface. During the first 10 4 to 10 6 years, the star cools mainly through neutrino emission. Later, 
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surface photon radiation dominates. The neutrino cooling rate of a neutron star is determined by 
the unknown state of matter in the core. Depending on the mass of the star and the properties 
of matter at several times nuclear density, nuclear matter may extend to the center of the star, 
or there may exist a distinct inner core consisting of a pion or kaon condensate, hyperon-rich 
matter, bulk quark matter, or some other exotic state. In standard cooling models, the neutrino 
emission is mainly through the modified Urea process taking place in the core. Later, neutrino 
bremsstrahlung in the inner crust may be important, though Pethick & Thorsson (1994) argue 
that band structure effects suppress this process. Neutrino bremsstrahlung processes in the core 
also contribute to the emission. In contrast with standard models, in accelerated cooling models 
neutrino emissivities are orders of magnitude larger. A number of plausible accelerated cooling 
mechanisms have been discussed, including, the direct Urea process on nucleons in a proton-rich 
core, free quarks, pions, or kaons. We summarize these neutrino processes in Table 1. 



2.2. Superfluid Effects 



Calculations of critical temperatures for nucleon pairing (see, e. g., Hoffberg et al. 1970; Chao, 
Clark, & Yang 1972; Ainsworth, Pines, & Wambach 1989; Wambach, Ainsworth, & Pines 1991) 
indicate that neutron stars form superfluids in their interiors shortly after their formation. In the 
inner crust, above the drip density of 4 x 10 11 g cm~ 3 , free neutrons form a 1 5o superfluid which 
flows through a lattice of neutron-rich nuclei. At or below nuclear density, p Q = 2.8 x 10 14 g cm~ 3 , 
the nuclei dissolve, forming a 3 P2 neutron superfluid, plus a dilute plasma of superconducting 
protons and normal electrons. (Recent work by Lorenz, Ravenhall, & Pethick [1993] and Lattimer 
& Swesty [1991] suggests nuclei dissolve at a lower density, 0.6p o .) Superfluid may or may not 
extend to the center of the star. 

Neutrino Emission and the Specific Heat: Superfluidity reduces the phase space volume 
available for excitations. Consequently, both the specific heat of the paired component, and any 
neutrino emission process that involves nucleons (see Table 1), would be suppressed. The direct 
Urea process on quarks would be unaffected, unless the quarks become superfluid. The possibility 
of quark pairing is unresolved; in this study we assume that quark superfluidity does not occur 
and there is no suppression of the quark Urea process. 

Frictional Heating: Friction between the normal matter of the star and a more rapidly 
rotating interior superfluid leads to heat generation which affects the star's thermal evolution. 
The strength and time-dependence of the heat generation depends on the microscopic physics of 
the interaction between the superfluid and the normal matter. We now turn to a discussion of this 
interaction, and derive the heating rate. 

The angular momentum of the rotating neutron superfluid resides in a velocity field containing 
singular regions of vorticity or vortex lines. The superfluid angular momentum is determined by 
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the density and distribution of these vortex lines, and can change only if they move relative to the 
star's rotation axis. Were the vortex lines immobilized, by the pinning process discussed below for 
example, the superfluid angular momentum would remain fixed. The rotational dynamics of the 
superfluid, and hence, the heating associated with the coupling between the superfluid and the 
crust, is determined by the motion of vortex lines. 

In early studies of neutron star thermal and rotational evolution, the core superfluid was 
assumed to be weakly coupled to the crust (Baym et al. 1969; Greenstein 1975; Harding, Guyer, 
& Greenstein 1978). However, subsequent investigations of neutron star spin fluctuations showed 
that at least 70% of the total moment of inertia of a neutron star is strongly coupled to the crust 
on time scales of greater than three days (Boynton & Deeter 1979; Boynton 1981). Theoretical 
studies suggest that the entire core is strongly coupled to the crust (Alpar, Langer, & Sauls 1984). 
Superconducting protons become entrained in the neutron circulation about a vortex, resulting 
in a supercurrent and an associated vortex magnetization. Electrons scatter efficiently off the 
vortex magnetic field, bringing the velocity between the core neutron superfluid and the charged 
component of the star into corotation over time scales of ~ 10 2 — 10 4 rotation periods (Alpar & 
Sauls 1988). On longer time scales, the core, which accounts for at least 90% of the total moment 
of inertia is tightly coupled to the crust. 

Unlike the core, the inner crust neutron superfluid, which makes up < 10% of the total 
moment of inertia, is only weakly coupled to the crust, allowing an angular velocity difference to 
develop between the superfluid and the crust as the pulsar spins down. Friction between the two 
components produces heat. The heating may occur steadily, if the superfluid and the crust slow 
at nearly the same rate, or suddenly, in a glitch. 

We now derive the internal heating rate in terms of the velocity difference between the 
superfluid and the rest of the star. As the star slows under an external torque iV cx t, the angular 
momentum of the inner crust superfluid, J s , changes according to 



where I c is the moment of inertia of the solid crust and the corotating core, J7 C is the spin rate of 
the crust, I is the total moment of inertia, and an overdot denotes a time derivative. The angular 
momentum of the inner crust superfluid in the Newtonian approximation is 



where r p is the distance from the rotation axis, p s is the superfluid mass density, and v s is the 
superfluid velocity, related to the angular velocity by v s (f,t) = r p £l s (r,t). The superfluid velocity 
changes as vortices move radially. Over evolutionary time scales, the superfluid and the crust 



i c n c + j s = N cxt = -i\a 
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remain close to rotational equilibrium such that Q s = Cl c = Q 



'OO- 



The total frictional heating rate is given by the rate of change of the total rotational energy 
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of the star, minus the rate at which work is done by the external torque: 

- |iVext|^ c = Pc\ J p s (r)r 2 p (Q s (r) - fl c )d 3 r, (3) 

where we have used eqs. (1) and (2) to obtain the final form. The integrand, which represents the 
heating rate per unit volume, is greatest near the equator and falls off near the rotational poles 
where r p vanishes. As we discuss below, the lag velocity between the superfluid and the crust, 
At; = r p ($l s — f2 c ), is nearly constant on spherical shells, except near the poles. Since we are 
studying the thermal evolution with a one-dimensional computer code, we average the superfluid 
heating on spherical shells to obtain 

Hit) = TT 2 \tl c \ J drp s (r)r 3 Av = \Q C \ J u cS dI = J eS \n c \, (4) 

where w e ff = 3-7rAt>/8r is the effective lag frequency and dl = (Ait /?>)p s (r)r A dr is the Newtonian 
approximation to the moment of inertia of a shell. This expression for the total heating rate is 
useful for comparison with the results of Shibizaki k, Lamb (1989) and USNT, who write the 
heating rate in this form. In our thermal evolution calculations, which allow for temperature 
gradients, we calculate the local heating rate in each zone, properly weighted by the moment of 
inertia with the general relativistic correction. 

The spin rate and its time-derivative are needed in the evaluation of the heating rate. We 
assume a magnetic dipole spin evolution, 

tl c (t) = -KS%, (5) 
where K is a positive constant. Integrating eq. (5) from a birth spin rate ri c (0) = f^o to Q c (t) gives 

n c (t) = "° • (6) 

yjl + 2KQ$t 

For times such that fi c (i) <C f^o, eq. (6) simplifies to 

Q c (t) = (2Kt)- 1 ' 2 = 190/i (jl-^j V2 s-\ (7) 

where we have chosen a fiducial value of the constant K that gives a spin down rate that is typical 
of the pulsars examined in this study. Each pulsar's spin evolution is characterized by a value of h 
(see Table 2). For the numerical solutions presented below, we use h = 1. For a birth period of ~ 
1 ms, eq. (7) is a good approximation after ~ 10 yr. 

The velocity lag appearing in the expression for the heating rate (eq. [4]) is determined by the 
vortex mobility, since only if vortices move radially can the superfluid velocity change. Interaction 
forces between vortices and nuclei, however, tend to pin vortices to the lattice. At high densities, 
above ~ 10 13 g cm -3 , the vortex-nucleus interaction is attractive, and vortices pin strongly to 
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nuclei (Alpar 1977; Epstein & Baym 1988). At lower densities, the interaction is repulsive, and 
vortices pin weakly to the interstices of the lattice (Epstein & Baym 1988; LE). In either region, 
vortex pinning tends to maintain a local superfmid rotational velocity that is higher than that of 
the crust. Following LEB, we consider the motion of vortices by the process of vortex creep, which 
proceeds by thermal or quantum activation of the vortices over their pinning barriers. Other 
mechanisms for moving vortex lines have been discussed. For example, Ruderman (1976; 1991) 
argues that in regions where vortex pinning is sufficiently strong, the crust may crack before the 
vortices unpin, and plates of crustal material will carry vortices toward the stellar equator. Jones 
(1990) has suggested that vortex pinning never occurs and that the superfluid remains in a state 
of near corotation with the crust. In LEB it was shown, however, that a vortex line adjusts its 
configuration in such a way as to be much more strongly pinned than Jones estimated. 

In this investigation, we assume that vortex creep allows the superfluid to slow at the same 
steady rate as the crust. Steady-state vortex creep, however, cannot account for all aspects of 
pulsar timing, especially the observed glitches in spin rate. Glitches are thought to reflect sudden 
transfers of angular momentum from the neutron superfluid to the crust. Some pulsars, e.g. the 
Vela pulsar, lose angular momentum mainly during glitches, rather than through steady vortex 
creep. In such cases, steady-state vortex creep is never achieved and the mean frictional heating 
rate is smaller than we calculate here. 

From vortex number conservation, the equation of motion for the superfluid is (see, e.g., LE) 

n s (f,t)^ ~ 2nsm v cr (f,t). (8) 
r p 

where v cr (r, t) is the radial vortex creep rate. In eq. (8) we have neglected gradients in Vl s , which 
are small in the inner crust. We consider cases such that \£l s — f2 c | <C £l c , and fix Cl s (r, t) = Cl c in 
cq. (8) to obtain 

v C r(r,t) ~ ~r P ^-- (9) 

The vortex creep rate has the general form 

v cr = v e- A ^, (10) 

where v is a microscopic velocity, A is the activation energy for a segment of vortex line to 
overcome its pinning barrier, and T c g- is the effective temperature and which differs from the stellar 
value at low temperatures because of quantum effects. The LEB result for the vortex creep rate 
is summarized in Appendix A. The activation energy depends sensitively on the vortex-nucleus 
interaction energy (which varies with the stellar density) and on the difference between the 
velocities of the superfluid and the crust, Av(r) = r p (Q s — fi c ). The solution of eqs. (9) and (10) 
for Av(r) depends only logarithmically on r p , so that Av{r) is nearly constant on spherical shells 
as we assumed in obtaining eq. (4). 

The strength of vortex pinning varies significantly with density in the inner crust, especially 
near 10 13 g cm~ 3 , where the character of pinning changes from interstitial to nuclear. An 
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interstitially vortex samples only the long range part of the interaction, and the pinning energies 
are very weak (~ 1 keV). Above ~ 10 13 g cm -3 , where vortices pin directly to nuclei, the pinning 
energies become difficult to calculate reliably, as the vortex samples the poorly understood 
short-range nuclear pinning potential. We consider two models for pinning which differ in the 
high-density region: 

• EB Pinning Model: LE and Epstein & Baym (1988) pinning parameters are used. 

• ACP Pinning Model: Above ~ 1.5 x 10 13 g cm -3 , pinning parameters are obtained in LEB 
from the approach of Alpar, Cheng & Pines (1989) based on the superfluid gap calculations 
of Ainsworth, Pines, & Wambach (1989). Elsewhere LE and Epstein & Baym (1988) pinning 
parameters are used. 

A calculation of the pinning interaction in the nuclear pinning region requires solving the 
BCS equations for pairing of neutrons in the presence of flow, both inside and outside nuclei. 
This problem has yet to be solved. In the EB model pinning energies are estimated by using 
Ginzburg-Landau theory; pinning energies in this model for the high-density region of the inner 
crust are typically ~ 10 MeV. In the ACP model, pinning energies are estimated by simple volume 
arguments; pinning energies in this model are ~ 1 MeV or less. 

Vortex creep will produce the most heat in the regions of nuclear pinning, since those regions 
have the most moment of inertia and the highest pinning energies. In the core, the velocity 
difference between components is negligible, and there is no significant frictional heat generation. 
Above ~ 1.3 x 10 14 g cm" 3 , the nuclear spacing and the characteristic vortex dimension become 
comparable. In this region the pinning may become "superweak", as suggested by Alpar et al. 
(1984), or may vanish altogether as Jones (1991) argues. In either case, the heating from such 
regions would be small and we ignore it in this investigation. 



3. SIMULATIONS 

We consider neutron stars with gravitational masses of IAMq constructed by solving the 
Oppenheimer-Volkoff equation for hydrostatic equilibrium based on the equations of state (EOSs) 
of Baym, Pethick, & Sutherland (1971; BPS; with the Pandharipande [1971] model C alternative 
at high density), Pandharipande & Smith (1975; PS), and Friedman & Pandharipande (1981; 
FP). The BPS model is representative of soft EOSs which yield stellar models with high central 
densities and thin, low-mass crusts. The PS model was chosen as an example of an extremely stiff 
model with low central density and a thick crust. The FP EOS, a representative intermediate 
model, results in a crust size within the range found using current models for the high density 
EOS (e. g., Prakash, Ainsworth, & Lattimer 1988). Radii R and central densities rho c for IAMq 
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stars constructed from these EOSs are given in Table 3. As a guide for comparison with Shibazaki 
& Lamb (1989) and UNST, we also show the maximum superfluid excess angular momentum, 
J = J dluj cr i t , and the maximum average lag, Co = J dl uj cvlt / J dl, integrated over the pinning 
region for our two pinning models. 

We consider two examples for the dominant core neutrino emission. In our standard cooling 
model, we use the modified Urea process, suppressed by superfluidity. As a specific example 
of accelerated neutrino emissivity, we consider the quark Urea process (Kiguchi & Sato, 1981). 
Electron-nucleus bremsstrahlung, the dominant neutrino emission process in the crust, is included 
in both cooling models. 

The frictional heating rate is obtained at each timestep as follows. Using the pinning 
parameters given as a function of density in Table 7, eq. (9), coupled with eqs. (7) and 
(Al)— (A15), is solved to obtain the superfluid velocity lag as a function of the local density and 
temperature. We evaluate the local heating rate |r2|cj c g<i/ on a grid of spherical mass shells, where 
dl is the moment of inertia of each shell, with the general relativistic correction. The thermal 
evolution is calculated using the penta-diagonal scheme of Van Riper (1991b), which allows for 
thermal gradients. Evolutionary calculations of this type (see also, e. g., Nomoto & Tsuruta 1987) 
show that isothermality is not obtained in the interior until an age of up to 10 4 yr, depending 
on the stellar model. Before the star becomes isothermal, the calculated luminosity can differ 
significantly from that found in an isothermal model. The neutrino emissivities, heat capacities, 
and thermal conductivities, and the details of the differencing scheme, solution algorithms, and 
boundary conditions are given in Van Riper (1991b). 

In Figs. 1-3 we show thermal evolution of neutron stars with and without frictional heating. 
The right ordinates give L^, the bolometric luminosity seen by a distant observer, and the 
left ordinates give the corresponding effective temperatures T s °°. Stars with thicker crusts, as 
illustrated by the PS model, have greater frictional heating and correspondingly higher surface 
temperatures at late times than stars with thinner crusts (e.g., BPS). Moreover, the EB pinning 
model, which leads to larger lag velocities and greater heating than the ACP model, produces 
larger surface temperatures. At late times the surface emission is controlled entirely by the heat 
source, and decreases as |fi c | oc i~ 3 / 2 [see eqs. (4) and (7)]. The temperature dependence of 
the vortex creep rate produces order unity deviations in the heating rate from the t~ 3 / 2 behavior, 
as shown in Fig. 4 (top panel). For comparison with the work of other authors, we present the 
scaled heating rate H(t)t 3 / 2 at 10 7 yr in Table 4 for different pinning models and equations of 
state. In Fig. 4 (bottom panel) we show the evolution of the inner crust temperature T m . The fall 
in T m after ~ 10 3 yr in the quark models produces a larger superfluid lag and a correspondingly 
larger heating rate. Comparison of our full evolutionary calculations with models that use 
time-independent energy generation are given in Appendix B. There we find that when significant 
frictional heating is present, the surface temperature is largely determined by balance between 
the source and the combined neutrino and photon emission, and is essentially independent of the 
star's thermal history. 
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3.1. Comparison with USNT 

USNT applied the heating formulae of Shibazaki & Lamb (1989) to thermal evolution models. 
Their heating is H(t) = J\Cl c \ where J is an assumed constant differential angular momentum. 
For each stellar model, three values of J are considered, corresponding to "superweak" , "weak" , 
and "strong" pinning, where «/ s trong = lO^weak = 100 ^superweak- Magnetic dipole breaking was 

_ 3 

assumed, giving H oc t 2 . 

The values of Superweak) ^strong' an< ^ the sca l e d heating rate H(t)t 3 / 2 are listed in Table 
5. USNT heating rates are plotted as a function of time in Fig. 5 for an FP star, together with 
our standard cooling model results. For PS stars, which have thicker crusts, our EB heating is 
~ 3 times stronger than shown in Fig. 5, and for a BPS star, which has a thinner crust, the EB 
heating is ~ 2 times less than shown. The changes for the ACP heating among different stellar 
models are relatively smaller. 

The results of USNT at early times (< 500 yr) differ from ours primarily because of the 
differing treatments of the early-time heating rate. USNT adopt a nearly constant heating rate of 
< 10 38 erg s _1 for times much earlier than 300 yr. However, superfluid condensation is expected to 
occur in the inner crust within a year of the star's birth, leading to substantial frictional heating at 
early times. Our treatment of early-time thermal evolution includes this time-dependent heating. 
The USNT cooling curves exhibit a small T at early times because their relatively low heating 
rate cannot compete with the rapid neutrino emission. The USNT results exhibit a rapid drop in 
surface temperature at ~ 30 yr caused by the arrival at the surface of the cooling wave from the 
core, rather than by the star's response to the changing heating rate. In our models, however, the 
temperature evolution has nearly zero slope at early times, when the heating and neutrino cooling 
rates essentially balance. The rapid drop in temperature does not occur until after ~ 100 yrs in 
our quark models because the initiation of the cooling wave is delayed by our larger early-time 
heating rate. At late times, when the thermal evolution is controlled by the same heating rate, 
our results are consistent with those of USNT. 



4. COMPARISON WITH OBSERVATIONS 

Soft x-ray emission has been detected from the surfaces of at least four neutron stars, and 
interesting upper limits have been obtained for numerous other objects. We now discuss the 
measurements, compare these data with our cooling simulations, and obtain constraints on the 
neutron star equation of state and internal heating source. The data are shown along with our 
cooling simulations in Figs. 1-3, and our results are summarized in Table 6. The interpretation of 
a spectrum in terms of an effective temperature T s °°, or equivalently, a luminosity L^, quantities 
we show in our figures, depends on the full spectrum of radiation emitted by the star. A first 
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estimate of the effective temperature could be obtained by assuming that the star radiates as a 
blackbody. However, the dependence of the atmospheric opacity on photon energy can lead to 
a substantial difference between the temperature inferred from a blackbody fit (Tbb) an d the 
effective temperature. For example, in a nonmagnetic atmosphere composed of light elements, 
the opacity drops with photon energy in the detector bandpass. As a result, the deeper, hotter 
atmospheric layers are exposed at higher energies, hardening the emitted spectrum, and the 
effective temperature can be up to a factor of three less than Tbb (Romani 1987). Strong 
magnetic fields broaden opacity features, softening the spectrum relative to the nonmagnetic 
case. Differences in treatment of the polarization-dependent transport have given different results. 
When considering only the portion of the spectrum in the ROSAT bandpass, Miller (1992; 1993) 
obtains effective temperatures about equal to Tbb for H, He, C, and N atmospheres, while for 
the EINSTEIN bandpass Tbb can be significantly less than the effective temperature — up 
to a factor of 0.7 in C and N atmospheres. On the other hand, Shibanov et al. (1993) obtain 
effective temperatures intermediate between Romani's (1987) nonmagnetic values and Tbb- in 
a non-magnetic Fe atmosphere, the average soft x-ray opacity is dominated by the presence of 
spectral lines, the spectrum is softened relative to a light-element atmosphere, and the effective 
temperature is close to Tbb (Romani 1987). The structure of heavy elements in strong magnetic 
fields is poorly understood, making an assessment of the full consequences for the spectrum 
impossible at present, however, because magnetic fields tend to soften the spectrum it is not 
impossible that such an atmosphere could have an effective temperature exceeding Tbb- Given the 
above caveats, in Figs. 1-3 we show the reported effective temperatures obtained from blackbody 
fits. 

Except for the Crab, whose age is known, we use the ages estimated from the pulsar spin 
down. This estimate is uncertain by at least a factor of three. For stars older than 100,000 
years, varying the age within this uncertainty may, depending on the model assumed, bring the 
observational point above or below our cooling curves. We show an error cross in the upper right 
of Figs. 1-3 to represent the uncertainties in the spindown ages and in the temperatures due to 
differences in atmosphere models. 

PSR 0656+14: The flux and spectrum over a limited band have been measured for this 
pulsar. From EINSTEIN observations, Cordova et al. (1989) deduce T s °° between 3 and 6 x 10 5 
K. Finley, Ogelman, & Kiziloglu (1992) obtain the spectrum with the ROSAT Position Sensitive 
Proportional Counter and find an effective surface black body temperature of T s °° ~ 9.0 ±0.4 x 10 5 
K. Radiation transfer in an atmosphere could have significant effects on the determination of 
T s °°. Finley, Ogelman, & Kiziloglu (1992) obtain ~ 2.4 ± 0.2 x 10 5 K using the non-magnetic He 
atmospheric model of Romani (1987). Finley & Ogelman (1993) find 9.2 ± 0.6 x 10 5 K using the 
magnetic (B = 9.4 x 10 11 G) He model of Miller (1992), while Anderson, et al. 1993, using the 
magnetic (B = 4.7 x 10 12 G ) H model of Shibanov et al. (1993), obtain 5.35 ± 0.35 x 10 5 K. 
Complicating the interpretation even further is a 14% pulsation in the thermal X-rays, suggesting 
a distribution of temperature over the surface. We show the blackbody results of Finley Ogelman, 



-12- 



k Kiziloglu (1992) in Figs. 1-3. 

Vela (PSR 0833-45): The flux and spectrum over a limited band have been measured for this 
star. The EXOSAT detections (Ogelman & Zimmermann 1989), fit assuming a blackbody source 
spectrum, are shown in in Figs. 1-3. 

PSR 1055-52: The flux and spectrum over a limited band have been measured by EXOSAT 
for this pulsar (Brinkmann & Ogelman 1987). We show the reported temperature error bar 
derived from a blackbody fit. 

Geminga: Halpern & Ruderman (1993) use a two-component model to fit the ROSAT 
spectrum, and obtain a high-temperature component of ~ 3 x 10 6 K, attributed to a polar cap, 
and a cooler component of ~ 5.2 ±1 x 10 5 K, associated with the remainder of the surface. The 
temperature of the cooler component is shown in Figs. 1-3. 

Crab (PSR 0531+21): We show an upper limit to for the Crab pulsar obtained from 
recent ROSAT measurements by Becker k Aschenbach (1994). This upper limit (and that for 
PSR 1929+10) is plotted at the same L^, rather than at the same T s °°, in Figs. 1-3, reflecting the 
experimental limit on the flux. 

PSR 1929+10: We show an upper limit to from this pulsar obtained by Alpar, et al. 
(1987). A recent measurement by Yancopoulos, Hamilton, k Helfand (1994) gives a comparable 
upper limit based on the net luminosity of the star. At the time and luminosity of this limit, 
the cooling curves are sensitive to the surface magnetic field strength chosen. Our curves assume 
B = 10 12 G. The curves would reach the PSR 1929+10 luminosity more quickly, by a factor of 
two, were B = 10 13 assumed; B = 10 11 would increase the cooling time by a factor of two. 

Other upper limits: The unlabeled short arrows in Fig. 1 are 2cr upper limits to the 
temperatures of radio pulsars from the ROSAT all sky survey (Becker, Triimper, k Ogelman, 
1993). These points were derived from observed counting rates assuming blackbody spectra. One 
ROSAT upper limit, for PSR 1916+14, is labeled in the figure because it provides interesting 
constraints on our models. 

Searches for neutron stars have proven unsuccessful in the following supernova remnants: 
Cas-A (Murray, et al. 1979), Kepler (Helfand, Chanan, k Novick 1980), Tycho (Gorenstein, 
Seward, k Tucker 1983), SN 1006 (Pye, et al. 1981), RCW 86, W28, (Helfand, Chanan, k Novick 
1980), and RCW 103 (Becker, et al. 1993; earlier interpretation of EINSTEIN data by Touhy k 
Garmire [1980] suggested a point source with a temperature T s °° = 2.2 x 10 6 K for RCW 103.) 
Present luminosity upper limits on point sources in Tycho, SN 1006, and RCW 103 fall below the 
standard cooling results; the discovery of a neutron star in any of these remnants would be strong 
evidence that accelerated neutrino cooling is taking place. 



5. DISCUSSION AND CONCLUSIONS 
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Frictional heating associated with differential rotation between the inner crust neutron 
superfluid and the crust could significantly affect neutron star thermal evolution, particularly at 
late times. Our simulations show that frictional heating keeps the surface temperature above 
10 5 K for t > 10 6 yr, perhaps to ages as great as 10 9 yr, suggesting a population of potentially 
observable stars (see Figs. 1-3). Stars with small radii and soft equations of state, as illustrated 
by our BPS model, are in closest agreement with all the current measurements and upper limits 
for standard cooling with ACP frictional heating or no heating (Figs. 1-3 nd Table 6). The 
stronger EB heating predicts a temperature at the age of PSR 1929+10 in excess of the upper 
limit for that pulsar. Our ACP model for a BPS star has an excess superfluid angular momentum 
of ~ 3.4 x 10 43 gm cm 2 rad s _1 , corresponding to an average lag of ~ 10 rad s~ 4 . Our conclusion 
differs from that of Alpar et al. (1987), who conclude that the excess angular momentum must 
be less than ~ 7 x 10 42 gm cm 2 rad s _1 , with a corresponding average lag of less than ~ 0.7 rad 
s~ 4 . Standard cooling of a BPS star lies slightly above the highest estimate for Vela's surface 
temperature. However, this discrepancy is not immitigable, given the uncertainties in Vela's age 
and emission spectrum. Because soft stars have thin crusts, frictional heating is less important 
than in stiff stars, and there is little difference between the cooling curves with EB heating, ACP 
heating, or no heating, except at very late times. Although soft equations of state are consistent 
with the neutron star surface emission, elsewhere we have argued that an equation of state as 
soft as BPS cannot account for pulsar post-glitch behavior (Link, Epstein, & Van Riper 1992). 
In addition, the maximum stable mass of a BPS star is only 1.414M Q , which is ~ 2% lower than 
the mass of PSR 1913+16. A star slightly stiffer than BPS, as preferred by Brown (1988), while 
consistent with both pulsar glitch data and neutron star mass determinations, would increase the 
nominal discrepancy between the cooling curves and the Vela detection. 

The Vela pulsar provides the strongest argument against standard cooling in stiffer stars, e. g. 
PS or FP, which have cooling curves lying above the Vela detection. Accelerated cooling models 
without heating fall below measured surface temperatures; the Vela and Geminga detections are 
the most compelling examples. Frictional heating lessens this disagreement, but the strongest 
frictional heating (e. g., EB in a PS star; Fig. 3) barely brings the cooling into agreement with 
the Vela and Geminga detections, while lying above the Crab and PSR 1929+10 upper limits. 
It appears that no combination of frictional heating and accelerated cooling is consistent with 
current data. Lattimer et al. (1993) and Page (1994), however, have emphasized that accelerated 
cooling processes are suppressed by superfluidity. Accounting for this suppression would give 
cooling curves that fall between those of our quark models and our standard models, and models 
of this sort may be consistent with current X-ray data. 
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APPENDIX A. VORTEX CREEP RATE 

The vortex mobility is determined by the strength with which vortices pin to the lattice, 
and by the local velocity difference Av = r p (O s — Q c ) between the superfluid and the crust. The 
pinning parameters used in this study are listed in Table 7 for selected densities: the pinning 
energy Uq, the vortex coherence length £, the Wigner-Seitz lattice spacing Z m i n , the normalized 
vortex stiffness r/A, and the critical differential velocity vb above which vortex lines cannot pin. 
LEB obtain the average creep velocity for vortex lines in two limits - the flexible limit (r <C 1), 
in which vortex motion occurs by unpinning at one site at a time, and the stiff limit (r » 1), in 
which a vortex segment unpins j bonds simultaneously. In these limits the creep rates are given by 

v i - ^imine" Al/Trff - 1 , flexible, (Al) 

1/2 

The symbols in eqs. (Al) and (A2) are defined as follows: 
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(A3) 
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(A5) 
(A6) 

.(j + l)WJ' (A7) 

where A = 1 — Av/vb and k = 0.0019 cm 2 s^ 1 is the quantum of vorticity. The quantities j and 
Aj are given by 

j = 3.242r 1 / 2 A- 1 / 4 (A8) 
Aj = 1.577l/oA 3/2 j, (A9) 
for A < 1/4, while for A > 1/4 they are 
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Here T = Av/vb = 1 — A, and, 
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The quantities r and j are functions of one another. In our calculations of the vortex 
creep steady state, these quantities are determined iteratively - a trial value of r is chosen, j is 
calculated, a new value of r is obtained, j is again calculated, etc., until convergence. 

Eqs. (Al) and (A2) are for the limiting cases of r <C 1 and r>l. To account for intermediate 
stiffness, we interpolate between the limiting cases by using the expression 




(A15) 
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APPENDIX B. COMPARISON WITH STEADY STATE SOURCE CALCULATIONS 

Here we compare the results of our evolutionary calculations with those of models that 
assume time-independent energy generation in the crust (e.g., Van Riper 1991a). We have 
calculated steady-state models for an inner crust source in an FP star. In Fig. 6 we plot the 
surface luminosity as a function of steady-state source strength. In Fig. 7 we compare the thermal 
evolution from our frictional heating model with the equilibrium temperature of a steady-state 
heating model. Each filled circle was obtained as follows: in our full evolutionary calculations we 
extract the heating rate at a selected time. We then perform an evolutionary calculation with 
this constant heating rate, until the temperature reaches a steady state. The circle represents the 
resulting surface luminosity. Where the circles fall on the cooling curves, heat generation balances 
the combined surface and neutrino emissions, and the thermal history of the star is unimportant. 

Cooling curves calculated assuming an isothermal interior are shown as dashed lines in Figs. 
6 and 7. These curves coincide with the evolutionary results for H < 10 36 erg s _1 for both 
standard and accelerated cooling, indicating that the interior temperature profile is effectively 
isothermal. The cooling curve for the standard cooling case is above the steady state points for 
10 35 > H > 10 32 erg s" 1 (10 4 to 10 6 yr). The additional radiation is residual thermal energy. For 
higher H, the hot crust insulates the interior from the surface and L is determined only by the 
heating rate in the crust (the initial heat content in the crust is negligible). The difference between 
the steady state and evolutionary models is greatest in our BPS model, which has relatively small 
heating, and is negligible for the highly heated PS star. The surface temperature of the accelerated 
cooling models is mainly due to the source once isothermality obtains. The heat content of the 
core, which has been reduced by enhanced neutrino emission, contributes only a negligible fraction 
to the surface luminosity. The quark model cooling curve shows a sudden decline between 100 and 
1000 years, similar to the drop caused by a cooling wave reaching the surface in a star undergoing 
accelerated cooling with no or weak crust heating. In our strong heating case, however, the 
shape of the cooling curve is controlled by the time dependence of the source rather than by the 
dynamics of conduction in the crust. 
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TABLE 1 

NEUTRINO EMISSION PROCESSES 







Emissivity a 




Process 


(erg s -1 cm -3 ) 


Reference 


Direct Urea -j 


f n — > p + e~ + v e 
I p + e~ — >■ n + z^ e 


~ 10 27 T 9 6 


Lattimer et al. 1991 


Modified Urea -j 


f n + n' — ► n' + p + e~ + P e 
[ n' + p + e~ — > n' + n + z/ e 
[ii + n'^ji + n' + i/ e + i/ e 


~ 10 20 T 9 8 


Priman & Maxwell 1979 


Core bremsstrahlung < 


. n + p^n+p + f e + P e 
[ e" + p — ► e~ + p + 1/ e + P e 


- 10 19 T 9 8 


Friman & Maxwell 1979 


Quark Urea -j 


f d — > it + e~ + z/ e 
I u + e~ — > d + z/ e 


~ 10 26 a c T 9 6 


Iwamoto 1980 


Pion condensate -j 


f n + 7r~ — ► n + e~ + v e 
I n + e~ — ► n + 7r~ + v e 


~ 10 26 T 9 6 


Maxwell et al. 1977 


Kaon condensate .... <j 


\ n + K~ ^ n + e~ + v e 
I n + e~ — ► n + + ^ e 


~ 10 24 T 9 6 


Brown et al. 1988 



Crust bremsstrahlung e + (A,Z) 

-► e - + (A, Z) + ^ e + ^ ~ 10 20 T 9 6 Itoh et al. 1984 



a Tg is the temperature in units of 10 9 K. 



TABLE 2 
SPIN DOWN PARAMETER 



Pulsar 


h 


Crab 


0.7 


Vela 


1.3 


Geminga 


2.6 


1055-52 


4 


1929+10 


8 


0656+14 


0.9 


1916+14 


0.3 
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Table 3 

M = 1.4M m NEUTRON STARS 





R 


Pc 


J/10 44 




EOS 


(km) 


(g cm -3 ) 


(g cm 2 rad s _1 ) 


(rad s _1 ) 


BPS 


7.34 


4.89 x 10 15 


4 a /0.34 6 


145 a /13 b 


FP 


10.57 


1.27 x 10 15 


24/2.4 


97/10 


PS 


15.83 


4.69 x 10 14 


450/19 


78/3.3 



a EB model 
ft ACP model 



TABLE 4 
HEATING RATES 

H(t)t 3 / 2 (10 Myr) 



EOS 


GAP 


(erg s 1 / 2 ) 


PS 


EB 


2.64 x 10 42 


PS 


ACP 


1.27 x 10 41 


FP 


EB 


1.50 x 10 41 


FP 


ACP 


1.52 x 10 40 


BPS 


EB 


2.55 x 10 40 



TABLE 5 
USNT HEATING RATES 



EOS 


J/10 44 


H(t)t 3 / 2 (10 Myr) 


Pinning^ 


(g cm 2 rad s _1 ) 


(erg s 1 / 2 ) 


PS/S 


31.0 


3.7 x 10 41 


PS/SW 


0.31 


3.7 x 10 39 


FP/S 


7.3 


6.4 x 10 40 


FP/SW 


0.073 


6.4 x 10 38 


BPS/S 


4.8 


2.9 x 10 40 


BPS/SW 


0.048 


2.9 x 10 38 



^S = strong pinning, SW = superweak pinning. 



TABLE 6 
ALLOWED CASES 



Model 



Pulsar 


Heating 


BPS 


BPS(accel) 


FP 


FP(accel) 


PS 


PS(accel) 


Prah 


Yl ("in 
IlUIlC 


v 

1 


Y 


Y 

1 


Y 


N 

1 N 


Y 




ATP 


1 


Y 


Y 


Y 


N 

IN 


Y 

1 




EB 


Y 


Y 


N 


Y 


N 


N 


Vpla 

V Lldi 




N 

1 N 


N 

in 


N 

IN 


N 

1 N 


N 

IN 


N 

IN 




ATP 


N 

1 N 


N 

IN 


N 

1 N 


N 

IN 


N 

1 N 


N 

IN 




EB 


N 


N 


N 


N 


N 


Y 


( . 1 ri era 

VJCIIIIII^Oj 




Y 


N 

IN 


Y 


N 

1 N 


Y 


N 

IN 




ATP 


Y 

1 


N 

IN 


N 

IN 


N 

IN 


N 

IN 


N 

IN 




EB 


Y 


N 


N 


N 


N 


Y 


1055-52 


none 


Y 


N 


Y 


N 


N 


N 




ACP 


Y 


N 


Y 


N 


Y 


N 




EB 


Y 


N 


Y 


N 


N 


N 


1929+10 


none 


Y 


Y 


Y 


Y 


Y 


Y 




ACP 


Y 


Y 


N 


N 


N 


N 




EB 


N 


N 


N 


N 


N 


N 


0656+14 


none 


Y 


N 


Y 


N 


Y 


N 




ACP 


Y 


N 


Y 


N 


Y 


N 




EB 


Y 


N 


N 


N 


N 


N 


1916+14 


none 


Y 


Y 


N 


Y 


Y 


Y 




ACP 


Y 


Y 


N 


Y 


N 


Y 




EB 


N 


Y 


N 


Y 


N 


Y 
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TABLE 7 

NUCLEAR PINNING PARAMETERS 



LOGio P 


U 






r/A 


V B 


( g cm~ 3 ) 


(MeV) 


(fm) 


(fm) 




(cm s _1 ) 




EB Pinnin 


g Model 




11.83 


0.00066 


2.6 


98. 


21. 


1.4 x 10 5 


11.99 


0.0018 


2.3 


93. 


20. 


1.6 x 10 5 


12.18 


0.0036 


2.0 


89. 


17. 


1.9 x 10 5 


12.41 


0.0090 


1.7 


84. 


16. 


2.1 x 10 5 


12.79 


0.025 


1.5 


79. 


13. 


2.7 x 10 5 


12.98 


0.4 


1.4 


71. 


0.053 


1.6 x 10 7 


13.18 


6.4 


1.5 


66. 


0.0068 


1.6 x 10 8 


13.53 


15.0 


1.9 


57. 


0.009 


3.8 x 10 7 


13.89 


9.0 


4.6 


39. 


0.070 


1.7 x 10 8 


14.12 


5.4 


18. 


29. 


0.12 


5.2 x 10 7 


ACP Pinning Model 


11.83 


0.00066 


2.6 


98. 


21. 


1.4 x 10 5 


11.99 


0.0018 


2.3 


93. 


20. 


1.6 x 10 5 


12.18 


0.0036 


2.0 


89. 


17. 


1.9 x 10 5 


12.41 


0.0090 


1.7 


84. 


16. 


2.1 x 10 5 


12.79 


0.025 


1.5 


79. 


13. 


2.7 x 10 5 


12.98 


0.4 


1.4 


71. 


0.053 


1.6 x 10 7 


13.18 


6.4 


1.5 


66. 


0.0068 


1.6 x 10 8 


13.53 


1.0 


7.9 


57. 


0.22 


8.0 x 10 6 


13.89 


0.97 


12. 


39. 


1.8 


3.2 x 10 6 


14.12 


0.22 


30. 


29. 


110. 


2.4 x 10 5 



- 21 - 



REFERENCES 

Ainsworth, T., Pines, D., k Wambach, J. 1989, Phys. Lett. B, 222, 173. 
Alpar, M. A. 1977, ApJ, 213, 527. 

Alpar, M. A., Anderson, P.W., Pines, D. k Shaham, J. 1984, ApJ, 276, 325. 
Alpar, M. A., Brinkmann, W., Kiziloglu, U., k Ogelman, H. 1987, A&A, 177, 101. 
Alpar, M. A., Cheng, K. S., k Pines, D. 1989, ApJ, 346, 823. 
Alpar, M. A., Langer, S. A., k Sauls, J. A. 1984, ApJ, 282, 533. 
Alpar, M. A. k Sauls, J. A. 1988, ApJ, 327, 723. 
Anderson, P. W., k Itoh, N. 1975, Nature, 256, 25. 

Anderson, S. B., Cordova, F. A., Pavlov, G. G., Robinson, C. R., k Thompson, Jr., R. J. 1993, 
ApJ, 414, 867. 

Baym, G., Pethick, C. J., Pines, D., k Ruderman, M. A. 1969, Nature, 224, 872. 
Baym, G., Pethick, C, k Sutherland, P. 1971, ApJ, 170, 299 (BPS). 

Becker, W. k Aschenbach, B. 1994, Lives of the Neutron Stars, Eds. A. Alpar, Kiziloglu, U., k 
van Paradijs, J. (Kluwer), in press. 

Becker, W., Triimper, J., Hasinger, G., k Aschenbach, B. 1993, in Isolated Pulsars, Eds. K. A. 
Van Riper, R. I. Epstein, k C. Ho (Cambridge: Cambridge), pll6. 

Becker, W., Triimper, J., k Ogelman, H. 1993, in Isolated Pulsars, Eds. K. A. Van Riper, R. I. 
Epstein, k C. Ho (Cambridge: Cambridge), pl04. 

Bildsten, L. k Epstein, R. I. 1989, ApJ, 342, 951. 

Boynton, P. E. 1981, in IAU Symposium 95, Pulsars, ed. R. Wielebinski k W. Sieber (Dordrecht: 
Reidel), p. 279. 

Boynton, P. E. k Deeter, J. E. 1979, in Compact Galactic X-Ray Sources, ed. F. K. Lamb k D. 
Pines (Urbana: University of Illinois), p. 168. 

Brinkmann, W. k Ogelman, H. 1987, A&A, 182, 71. 

Brown, G. E. 1988, Nature, 336, 519. 

Brown, G. E., Kubordera, K., Page, D., k Pizzochero, P. 1988, Phys.Rev.D, 37, 2042. 
Chao, N.-C, Clark, J. W., k Yang, C.-H. 1972, Nuc.Phys.A, 179, 320. 
Chau, H. F. k Cheng, K. S. 1993a, Phys.Rev.B, 47, 2707 

Chau, H. F. k Cheng, K. S. 1993b, in Isolated Pulsars, Eds. K. A. Van Riper, R. I. Epstein, k C. 
Ho (Cambridge: Cambridge), p35. 

Cordova, F. A., Hjellming, R. M., Mason, K. O., k Middleditch, J. 1989, ApJ, 345, 451. 



- 22 - 



Epstein, R. I., & Baym, G. 1988, ApJ, 328, 680. 

Finley, J. P. k Ogelman, H. 1993, in Isolated Pulsars, Eds. K. A. Van Riper, R. I. Epstein, k C. 
Ho (Cambridge: Cambridge), pllO. 

Finley, J. P., Ogelman, H., k Kiziloglu, U. 1992, ApJ, 394, L21. 

Friedman, B., k Pandharipande, V. R. 1981, Nuc.Phys.A, 361, 502 (FP). 

Friman, B. L. k Maxwell, O. V. 1979, ApJ, 232, 541. 

Gorenstein, P., Seward, F., k Tucker, W. 1983, Helfand, D. J. 1983, in IAU Symposium 101, 
Supernova Remnants and Their X-Ray Emission, Ed. P. Gorenstein k J. Danziger (Dordrecht: 
Reidel), p. 1. 

Greenstein, G. 1975, ApJ, 200, 281. 

Halpern, J. P., k Ruderman, M. 1993, ApJ, 415, 286. 

Harding, D. Guyer, R. A., k Greenstein, G. 1978, ApJ, 222, 991. 

Helfand, D. J. 1983, in IAU Symposium 101, Supernova Remnants and Their X-Ray Emission, 
Ed. P. Gorenstein k J. Danziger (Dordrecht: Reidel), p. 471. 

Helfand, D. J., Chanan, G. A., k Novick, R. 1980, Nature, 283, 337. 

Hoffberg, M., Glassgold, A. E., Richardson, R. W., k Ruderman, M. 1970, Phys.Rev.Lett, 24, 775. 
Itoh, N., Kohyama, J., Matsumoto, N., k Seki, M. 1984, ApJ, 285, 304. 
Iwamoto, N. 1980, Phys.Rev.Lett, 44, 1637. 
Jones, P. B. 1990, MNRAS, 243, 257. 
Jones, P. B. 1991, ApJ, 373, 208. 

Kiguchi, M., k Sato, K. 1981, Progr.Theor.Phys.Lett, 66, 725. 

Lattimer, J. M., Pethick, C. J., Prakash, M., k Haensel, P. 1991, Phys.Rev.Lett, 66, 2701. 
Lattimer, J. M. k Swesty, F. D. 1991, Nuc.Phys.A, 535, 331. 

Lattimer, J. M. k Van Riper, K. A., Prakash, M., k Prakash, M. 1993, ApJ425, 802. 

Link, B. k Epstein, R. I. 1991, ApJ, 373, 592 (LE). 

Link, B., Epstein, R. I., k Baym, G. 1993, ApJ, 403, 285 (LEB). 

Link, B., Epstein, R. I., k Van Riper, K. A. 1992, Nature, 359, 616. 

Lorenz, Ravenhall, D. G., k Pethick, C. J. 1993, Phys.Rev.Lett, 70, 379 (1993) 

Maxwell, O. V., Brown, G. E., Campbell, D. K., Dashen, R. F., k Manassah, J. T. 1977, ApJ, 
216, 77. 

Migdal, A. B. 1959, Nuc.Phys, 13, 655. 
Miller, M. C. 1992, MNRAS, 255, 129. 



-23- 



Miller, M. C. 1993, in Isolated Pulsars, Eds. K. A. Van Riper, R. I. Epstein, k C. Ho (Cambridge: 
Cambridge), pl53. 

Murray, S. S., Fabbiano, C, Fabian, A. C, Epstein, A., k Giacconi, R. 1979, ApJ, 234, L69. 

Negele, J. W. k Vautherin, D. 1973, Nuc.Phys.A, 207, 298. 

Nomoto, K. k Tsuruta, S. 1987, ApJ, 312, 711. 

Ogelman, H. k Zimmermann, H.-U. 1989, A&A, 214, 179. 

Page, D. 1994, ApJ, 428, 250. 

Pandharipande, V. 1971, Nuc.Phys.A., 178, 123. 

Pandharipande, V. R., k Smith, R. A. 1975, Nuc.Phys.A, 237, 507 (PS). 
Pethick, C. J. k Thorsson, V. 1994, Phys. Rev. Lett, in press. 

Pines, D., Shaham, J., Alpar, M. A., k Anderson, P. W. 1980, Prog.Theor.Phys.Suppl., 69, 376. 

Prakash, M., Ainsworth, T. L., k Lattimer, J. M. 1988, Phys.Rev.Lett, 61, 2518. 

Pye, J. P., Pounds, K. A., Rolf, D. P., Seward, F. D., Smith, A., k Willingale, R. 1981, MNRAS, 
194, 569. 

Romani, R. W. 1987, ApJ, 313, 718. 
Ruderman, M. 1976, ApJ, 203, 213. 
Ruderman, M. 1991, ApJ, 366, 261. 

Shibanov, Yu. A., Zavlin V. E., Pavlov, G. G., Ventura, J., k Potekhin, A. Yu. 1993, in Isolated 
Pulsars, Eds. K. A. Van Riper, R. I. Epstein, k C. Ho (Cambridge: Cambridge), pl74. 

Shibazaki, N. k Lamb, F. K. 1989, ApJ, 346, 808. 

Tuohy, L, k Garmire, G. 1980, ApJ, 239, L107 

Van Riper, K. A. 1988, ApJ, 329, 339. 

-1991a, ApJ, 372, 251. 

-1991b, ApJS, 75, 449. 

Van Riper, K. A., Epstein, R. l.,k Miller, G. S. 1991, ApJ, 381, L47. 

Umeda, H., Shibazaki, N., Nomoto, K., k Tsuruta, S. 1993, ApJ, 408, 186 (USNT). 

Wambach, J., Ainsworth, T. L., k Pines, D. 1991, in Neutron Stars: Theory and Observations, 
ed. J. Ventura k D. Pines (Dordrecht: Kluwer), 37. 

Yancopoulos, S., Hamilton, T. T., k Helfand, D. J. 1994, ApJ, (in press). 



This preprint was prepared with the AAS IATgX macros v3.0. 



-24- 



Fig. 1. — Thermal evolution of an FP star of gravitational mass IAMq without frictional heating 
(solid curves) and with heating using the EB (dashed curves) and ACP (dot-dashed curves) pinning 
parameters. For each case, the evolution was calculated with "standard" neutrino emissivities 
and with "accelerated" cooling due to emission from quarks. The dotted curve is an isothermal 
approximation to the standard EB model. A surface magnetic field of 10 12 G is assumed. Data 
and upper limits, discussed in the text, are shown. Uncertainties associated with the atmospheric 
physics and age determinations are shown by the error cross on the right side of each figure (see 
text). When comparing theoretical curves with data, two factors must be borne in mind. 1) For all 
pulsars except the Crab, the ages shown correspond to a breaking index of 3; an index of 2 gives 
half the age shown, while an index of 4 give 1.5 times the age shown. 2) The theoretical curves 
correspond to an h = 1 spin-down parameter (see eq. [7]) The heating rate of a given pulsar is 
proportional to its spin-down parameter (given in Table 2). At late times, the luminosity equals 
the heating rate, and the curves with heating must be shifted by a factor h in the luminosity. At 
early times, the required shift can be estimated by the difference between the EB and ACP curves, 
which differ by a factor of ~ 10 in the heating rate. 

Fig. 2. — Same as Fig. 1, for a PS star. 

Fig. 3. — Same as Fig. 1 for a BPS star, except accelerated cooling with heating was calculated 
assuming an isothermal core to save computer time. 

Fig. 4. — Evolution of the scaled heating rate, H(t)t 3 / 2 , normalized to its value at at 10 7 yr (top 
panel), and the inner crust temperature (bottom panel). 

Fig. 5. — Net heating rates for an FP star from this work (solid lines) and from USNT (dashed 
lines) for several pinning models. 

Fig. 6. — Steady state surface luminosity as a function of inner crust source rate for FP neutron 
stars. The luminosity is equal to the source rate along the dotted line. 

Fig. 7. — Cooling curves for an FP star compared with steady state results (filled circles). The 
logarithm of the source strengths H in erg s _1 is shown above the bottom axis. Models calculated 
assuming isothermality are shown as dashed lines. 
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